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We study both experimentally and numerically the steady zonal flow generated by longi- 
tudinal librations of a spher i cal ro tating container. This study follows the recent weakly 
nonlinear analysis of Bussel ( 2010t ). developed in the limit of small libration frequency 



- rotation rate ratio, and large libration frequency - spin-up time product. Using PIV 
measurements as well as results from axisymmetric numerical simulations, we confirm 
quantitatively the main features of Busse's analytical solution: the zonal flow takes the 
form of a retrograde solid body rotation in the fluid interior, which does not depend on 
the libration frequency nor on the Ekman number, and which varies as the square of the 
amplitude of excitation. We also report the presence of an unpredicted prograde flow at 
the equator near the outer wall. 



1. Introduction 

Longitudinal librations (reported below as librations) are periodic oscillations of a 
rotating container about its axis of rotation. Despite the fact that theses oscillations 
are time-dependent, it has been recently suggested that they can generate non-linear ly 
a steady axi symmetric flow in the liquid interior through the Ekman boundary layer 
( Bussel l20ld) . A better know ledge of this resu lting flow is of great interest in geo- and 



astrophysics (see for instance Noir et al. 20091) where libration is driven by gravitational 



interactions and is used to in vestigate the interior structure of planets (e.g. iMargot et al 
l2007tlVan Hoolst et alfooOht 



Despite of the possible applic ations, flows driven by lib ration in rotating contain- 
ers have not been much studied. Aldridge fc Toomre ( 19691) have observed experimen- 



tally that inertial modes can be excited by li bration at parti cular resonance frequen- 
cies, which has been confir med numerically by iRieutordl (|1991l ) . However in the case of 
Aldridge &: Toomrd ( 19691 ). the experimental results are measurements of pressure dif- 



ferences between two points on the axis of rotation , and d o not bring information about 
the resulting flow created in the interior. iTilgnerl (|1999l) has investigated numerically 
the linear response to the forcing in the case of a spherical shell and has shown that 
the pres ence of an i nner c ore only marginally modifies the resonance frequencies. More 
recently iNoir et al\ (|2009h have studied experimentally by direct flow visualization the 
presence of centrifugal instabilities in the form of Taylor-Gortler vortices near the outer 
boundary, by varying the frequency and the amplitude of libration. The same group has 
also performed LDV measurements of libration-driven zonal flows in a librating cylinder 
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( Noir et al. 2010| ) i n the case of high fr equency librations and axisymmetric simulations 
in a spherical shell ( Calkins et al. 2010). Finally, a complete weakly non-linear theory of 
th e zonal flow d riven by low frequency librations in a sphere has been recently developed 
bv iBussd <l2010h in the absence of direct resonant forcing of any inertial wave. To the 
best of our knowledge, the main features of this analytical solution have not yet been 
validated quantitatively. This is the aim of the present work, combining experimental 
and numerical approaches. The paper is organized as follows. §2 gi ves a b r ief su mmary 
of the governing equations and of the weakly nonlinear analysis of Bussd (|201dh . In §3 
we present the experimental setup and the numerical model used in this study. Then 
experimental and numerical results are compared to the theory in §4. Discussion and 
conclusion are given in §5. 



2. Weakly nonlinear theory 

Let us consider a spherical cavity of radius R filled with a homogeneous and incom- 
pressible fluid of kinematic viscosity v. In the inertial frame the cavity rotates with an 
angular velocity 

n(t) = ( fi + cos(wia t)) k, (2.1) 
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where fio is the mean rotation rate, Af2 is the amplitude of libration, wub the libration 
frequency and k the unit vector in the direction of the rotation axis. Using respectively 
R and f2o _1 &s lengthscale and timescale, the dimensionless equations of motion written 
in the frame rotating at the angular velocity Qq and the sidewall boundary conditions 
are given by 

^ + u ■ Vit + 2k x u = -Vp + EV 2 u (2.2a) 
ot 

V • u = (2.26) 

u = e k x r cos(w t) at \r\ = 1, (2.2c) 

where r is the spherical radial coordinate, u the velocity measured in the rotating frame, 
p the modified pressure taking into account centrifugal forces, E = v/SIqR 2 the Ekman 
number, e = A£!/2f2o the normalized amplitude of libration and w = uJub/Qo the nor- 
malized librational frequency. The librational forcing appears in the problem throu gh the 
boundary condition (12.2b 1 ). This problem has been recently solved bv lBussel ( 201dh in the 
limit: 

\/£«u;<6<1. (2.3) 

The limit of small Ekman number allows splitting of the velocity field into two parts: a 
component U describing the inviscid flow in the interior and a boundary layer component 
u. Following the weakly non-li near method previously us ed in the case of precession 
(|Busselll968l i and tidal forcing (jSuessI Il97ll ). IBussd (|2010h finds an expression for the 



steady zonal flow in the inviscid interior U in the limit of low libration frequencies 

U = e 2 k x rf(\k x r\ 2 ) (2.4a) 
„ 2N 259a; 2 - 360 

where f ^ = 2400(1 -*y (2 ' 46) 

The function f(x 2 ) represents the average difference in angular velocities between the 
container and the fluid divided by e 2 . For < r < 1, f{x 2 ) is negative, i.e. the fluid 




Figure 1. (a) Photo and (b) sketch of the experimental setup. 

is expected to rotate in the retrograde direction. Moreover this differential rotation is 
nearly constant up to r ~ 0.6 with a mean value of —0.154. The zonal flow can thus be 
assimilated in the bulk to a retrograde solid body rotation superimposed on the mean 
rotation, whose amplitude is independent of the libration frequency and Ekman number 
and changes as — 0.154 e 2 . In the present paper, we verify experimentally and numerically 
these main features. Note that (j2.41 b) diverges for x = 1, i.e. near the outer boundary at 
the equator. Here, the analytical approach requires t he introduct ion of a specific scaling 
due to the singularity of the Ekman boundary layer ( Bussei 2010h . 



3. Methods 

3.1. Experimental setup 

Figure [T] shows a photo and a schematic vie w of the experimen tal setup used in this 
study, which is the same as the one used by Morize et al. ( 2010T ) to study zonal flows 
driven by tides. It consists in a hollow sphere, of radius R = 10 cm, which was molded 
in a transparent silicone gel to allow flow visualization. The sphere is filled with water 
and seeded with Optimage particles of 100 fim in diameter and of density 1 g.cm~ 3 ±2%. 
The sphere is set in rotation about its vertical axis {Oz) with a mean angular velocity 
£lo up to 85 r.p.m with a precision of ±0.3%. Once a solid body rotation is reached 
(typically in ~ 10 minutes), a librational motion is set using sinusoidal oscillations of the 
angular velocity of the sphere of the form e cos(w;ib t) where ujm, can be chosen between 
0.6 and 120 r.p.m with a precision of ±0.3%. In terms of dimensionless numbers, we have 
explored the following ranges: the Ekman number E = v/SIqR 2 £ [10 -5 ; 10~ 4 ], the ratio 
between the libration and spin frequency u> — wjib/^o £E [0.04; 0.1] and the amplitude of 
libration e £[0.02; 0.15] with a precision of ±0.6%. 

In order to measure the velocity field in the equatorial plane induced by librational 
forcing, we used a rotating particle image velocimetry (PIV) system. A miniature wire- 
less camera 1/4' Sharp HighQ CCD 29.4 x 22 mm of resolution 576 x 768 pixels rotates 
at a constant angular velocity Qq and measurements are made from above through the 
transparent top surface. The PIV particles are illuminated by a laser sheet, of thick- 
ness 3 mm, produced by a continuous laser (4 W) in the equatorial plane. We wait for 
about 20 oscillations after turning on the libration forcing to ensure that the response 
of the fluid is well established, then we start acquiring pictures for PIV measurements 
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using a video transmitter -receptor system. Velocity fields are computed using DPIVSoft 
(jMeunier fc Leweke 20031 ) on a 60 x 80 grid with a spatial resolution of 3 mm, close to 



the laser sheet thickness. We look for the time independent axisymmetric zonal flow in- 
duced by the libration whereas the forcing of the sphere is of the form e cos(a;t). Velocity 
fields are thus time-averaged over several periods of libration in order to eliminate the 
time dependent term. This also significantly enhances the signal to noise ratio. However, 
our experimental setup only allows collection of PIV data for a limited time. Hence, we 
cannot set the libration frequency at too low a value because we would not be able to 
average out the resulting data over enough periods to have a correct velocity profile. In 
the experiments we have considered frequencies in the range 0.04 u) ^ 0.1. Higher and 
lower libration frequencies have been studied using numerical methods, that we describe 
in the following section. 

3.2. Numerical approach 

In addition to the experiments, we have performed axisymmetric numerical simulations 
of the flow within a sphere of radius R in rotation with an angular velocity f2(£) = 
£lo(l + e cos(cj£))e z . We use a commercial software, Comsol Multyphisics©, based on 
the finite elements method to solve this problem. The numerical grid consists of two 
domains: (i) a boundary layer domain of thickness 0.035 R all along the outer boundary 
and the axisymmetric axis, which is discretiscd in the direction normal to the boundary 
into 25 quadrilateral elements with initial thickness of 5.10 _5 i? and a stretching factor 
of 1.2; (ii) a bulk zone with triangular elements. All elements are of standard Lagrange 
PI — P2 type (i.e. linear for the pressure field and quadratic for the velocity field). Note 
that the finite element method does not induce any particular problem around r = 
and t hat no stabilization tech nique has been used in this work. The temporal solver is 



IDA (jHindmarsh et all 12005.1 ). based on backward differencing formulas. At each time 
step the system is solved with the sparse direct linear solver PARDISC0 The number 
of degrees of freedom (DoF) used in the simulations is constant and equal to 157 869 
DoF. Our numerical model solves the Navier-Stokes equations in the frame rotating at 
the velocity £!o e z, with no-slip boundary conditions and a fluid initially at rest in this 
frame (i.e. a solid body rotation at fio in the inertial frame). At time t — 0, libration of 
the outer boundary is turned on and computations are pursued until a stationary state 
is obtained, which is reached typically in less than 10 libration periods. The velocity is 
then averaged in time over 5 libration periods to obtain the steady zonal flow. Results are 
non-dimensionalised as in the experiment and in the the ory. The numerical model has 
been validated in reproducing the experimental results of Aldridge fc Toomre ( 19691) . In 



their paper they define a fixed libration Reynolds number Re^ = uJubR 2 jv = 6.2 ■ 10 4 
and their applied angular velocity is given by 

O(t) = Sl + e uub cos(uJubt) where e = 8.0 7r/180 rad. (3.1) 

Pressure measurements from our numerical simulation when systematically changing 
the libration frequency w are pr esented in figure [21 a and sh ow an excellent agreement 



wi th the ex p erime ntal results of Aldridge fc Toomre ( 19691 ) and the numerical results 



of Rieutordl (1991), which validates the numerical model. In figure EJa, each peak cor- 
responds to the resonant forcing of a given inertial mode o f the rotating sphere and is 
labeled by two integers (n,m)(see Aldridge fc Toomrel 1969f ). There is no inertial mode 



for | a; | ^ 2 and the modes are progressively damped when 1/w increases due to the re- 
duced coupling between the container's oscillation and the fluid interior as well as to the 
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Figure 2. (a) Time-averaged crest-to-trough amplitude of pressure difference C p between the 
center and the pole of the sphere for various frequency ratios w = uJub/^lo and for a libration 
forcing ^lo + euub cos(wj<&t) with I = 8.0 ir/l80rad and Re^ — LJubF t 1 Jv = 6.2- 10 4 . The squares 



stand for our numerical values, the circles for experimental values of lAldridge fc Toomrd (|1969f ) 
and the line for the theoretical plot, (b) Velocity in the z-direction (i.e. in the direction of the 
rotati on axis) at time uit = 3tt/2 [2n] for the mode (2, 1) in the notation of Aldridgc fc Toomrd 
(1969), corresponding to l/u = 1.066. 
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Figure 3. (a) Time-averaged velocity field obtained by PIV measurement in the equatorial 
plane for E — 1.15 x 10 -5 , e = 0.08 and lo — 0.1. The background is colored as the norm of the 
horizontal velocity. The center of the sphere is at (0,0). (b) Mean experimental dimensionless 
azimuthal velocity (sym bols) correspo nding to the velocity field of Fig. [3] a and comparison with 
the theoretical results of iBussel (|2010f ) (dashed-dotted line). 



increased viscous damping as the structure of the forced modes bec omes more complex. 
The v elocity in the z-direction for the mode (2, 1) in the notation of lAldridge fc Toornre 
(Il969l) is presented in figure [2jb and shows the inertial wave excited by libration forcing 
as wel l as the structure of the outer boundary layer. In the next section, following iBusse 
(201fj), we investigate the limit lu <C 1 where the forcing of inertial modes is negligible, 
but where a global zonal flow is excited. 



4. Results 

Figure [3ja shows an example of the velocity field obtained by PIV measurement in 
the equatorial plane. The stationary flow is azimuthal and axisymmetric. Besides, as the 
system rotation is clockwise, the zonal flow corresponds here to a retrograde circulation 
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Figure 4. Dimensionless time-averaged velocity profiles obtained by numerical simulation in 
the equatorial plane with E — 5 x 1CP 5 , e = 0.2 and u) = 0.03 (dashed line), uj = 0.06 ( dotted 
line) a nd ui = 0.1 (dashed-dotted line), compared with the inviscid analytical solution of[Bussc 
(2010) (continuous line). 



opposed to the rotation of the sphere. From this velocity field we can plot an averaged 
azimuthal velocity profile at the equator for given parameters in term of dimensionless 
quantities as a function of the radial distance (Fig. (3]b). Experimentally, due to optical 
deformation induced by the planar air-silicone and spheroidal silicone-water interfaces, 
it is not possible to measure the profile for r > 0.85. We observe the steady zonal flow 
already visible in figure|3]a. Moreove r we can direc tly compare this dimensionless quantity 
with the analytical solution given bv lBussd <l2010h and we observe an excellent agreement 
up to r ~ 0.6 with no adjustment parameter. An example of the velocity field obtained 
numerically is shown in figure|4l which also exhibits good agreement with the analytical 
solution in the bulk. For r ^ 0.6 a deviation in the prograde direction with respect to the 
theoritical profile due to the librating outer boundary is observed and will be discussed 
below. But for now, we concentrate on the mean zonal flow induced in the bulk. 

We have performed series of experiments and numerical calculations to systematically 
check the effect of the three control parameters E, e and uj on this bulk zonal flow. To 
do so, we define a reproducible method to synthetise the experimental and numerical 
data. Since the function f{x 2 ) (|2.4l b) may be considered as constant up to r ~ 0.6, 
which means that the predicted zonal flow is almost a solid body rotation up to r ~ 0.6, 
we take the average value of the measured \U\/r between r = 0.1 and r — 0.6, i.e. the 
non-dimensionnalized mean angular velocity, and compare it with the theoretical value 
0.154e 2 . Note that experimental results are represented by bars (see for example figure 
[5]) which represent both the uncertainties of the PIV measurements and the deviation of 
the measured velocity profile from a pure solid body rotation. 

In figure [5] we investigate the influence of the amplitude of libration e on the zonal 
flow. Experimentally we set E = 2.3 x 10~ 5 , uj = 0.07 and we systematically change e 
between 0.02 and 0.15. To explore a larger range of amplitude we have also performed 
numerical simulations with E = 4 x 10~ 5 , uj = 0.04 and e €[0.01; 0.2]. Both experimental 
and numerical results are quantitatively compatible with the theory with no adjustment 
parameter. The steady azimuthal velocity scales as e 2 for a large range of e. When e 
becomes larger than 0.2 the weakly nonlinear hypothesis cannot be used anymore because 
terms of higher order cannot be neglected. We also notice that even in a range of values 
where the condition uj <C e in (|2.3[) is not fully satisfied, the zonal wind intensity still 
scales as e 2 . In fact, rather than the more restrictive condition (|2.3|) . we only require that 
(i) E <C 1 to decouple the bulk and boundary layer flows, (ii) \f~E < w < 1 in order 
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Figure 5. Dimensionless average amplitude of \U\/r between r — 0.1 and r — 0.6 as a function 
of the amplitude e. Experimental results (bars) for E = 2.3 x 10~ J , uj = 0.07, and nu merical 
result s (squares) for E — 4 x 10 -5 , uj = 0.04 are compared with the theoretical value of|Bussc 
l|2010h (dashed-dotted line). 
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Figure 6. Dimensionless average amplitude of \U\/re 2 between r — 0.1 and r = 0.6 as a function 
of the Ekman number. Bars are experimental results for e — 0.08 and uj = 0.1 and squares are 
n umerica l resul ts for e = 0.2 and oj = 0.06. The dashed-dotted line shows the theoretical result 
of iBussel (2010). The velocity profiles have been rescaled by e 2 following the results presented 
in figure [S] Experimentally, errorbars increase with the Ekman number mainly due to the fact 
that our experimental setup does not allow us to average the velocity on a sufficient number of 
periods when increasing E. 



to neglect the excitation of inertial waves ( Aldridge k, Toomre 1969) and to ensure that 
the spin-up effect of the libration is confined inside the outer boundary layer, and (iii) 
e< 1 to remain in the weakly nonlinear regime. 

In figure[6]we report the systematic study of the influence of the Ekman number on the 
zo nal flow. Expe rimental results are compatible with the no-Ekman dependence predicted 
bv lBussd (|2010l) in a large range of Ekman number with no adjustment parameter. This 
is confirmed numerically up to E ~ 10~ 3 . For larger values of the Ekman number, 
the condition (|2.3I) is not fulfilled and we cannot assume that the effect of spin-up is 
negligible in the bulk. We have also noticed numerically that further decreasing E for a 
given to = 0.1 (which is not so small) leads to a deviation from the theory. Indeed, forced 
inertial modes are not expected to be negligible anymore and can perturb the zonal flow. 
In particul ar, nonlinear self-in teraction of these forced mode can drive localised zonal 
winds (e.g. iMorize et al. 2010l ). These peculiar behavior appears out of the asymptotic 
limit (|2.3p under consideration here and will be the subject of a future study. 

Figure [TJa shows that the amplitude of the zonal wind in the bulk does not depend on 
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Figure 7. (a) Dimensionless average amplitude of \U\/re 2 between r — 0.1 and r — 0.6 as a 
function of the frequency of libration uj; squares are numerical values (E = 5 X 10~ 5 , e = 0.2) 
and bars are experi mental data (E — 1.5 x 10~ 5 , e = 0.1). The dashed-dotted line shows the 
theoretical result of Bussc (2010). The velocity profiles have been rescaled by e 2 following the 
results presented in figure [5] (b) Distance of the minimum of the velocity profile from the outer 
boundary as a function of ui obtained by numerical simulation for E = 5 x 10~ 5 , e = 0.2. The 
dotted line scales as 1/y/u, representative of a skin effect. 



oj, as suggested by Busse. Nevertheless, as can be seen in figurelU the flow near the outer 
wall changes with oj. We expect this prograde flow to be related to the same mechanism 
of bounda r y laye r ejection near the critical latitude as the prograde jets described in 
iNoir et all (2009). To quantify the distance at which the real flow deviates from the 
analytical solution, we identify the value r m i n where the velocity has a minimum. The 
thickness of the layer where the prograde flow develops is plotted in FigJTlb and is found 
to scales as l/y/w, which is representative of a skin effect. So if oj becomes too small at a 
fixed E, the layer where the effects of external walls are important is visible in the bulk 
and perturbs the zonal flow. 



5. Conclusion 

In this paper, combining numerical and experimental studies, we report the first quan- 
titative measurements of the steady flow driven by longitudinal librations in a ro tating 
sphere . This approach confirms the main features of the weakly nonlinear theory of lBusse 
( 20101 ): a retrograde differential rotation induced by the libration of the sphere takes place, 
which may be assimilated to a solid body rotation for r < 0.6. It is also shown that the 
amplitude of this steady zonal flow is independent of w and E and scales as e 2 . Note that 



the s ame features have been observed experimenta lly in a librating cyl inder (jNoir et al 



20101) and numerically in a librating spherical shell (jCalkins et aZ.|[2010h and thus appear 

to be generic of librating flows. 

The main differences between our results and the theoretical profile of iBussd (|2010l) 
arise close to the outer boundary at the equator. There, we observe a prograde flow 
in a layer of thickness proportional to 1/y/oj. The analytical resolution of this peculiar 
feature would necessitate a special trea tment since it appears at the critical latitude of 
the outer boundary layer (Busse |2010| ). This, as well as the experimental study of a 
librating spherical shell, will be the subject of forthcoming studies. 
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